3-D Laplace Equation And Spherical Harmonic Function

1 3D Laplace Equation

三维拉普拉斯方程方程可记成多种形式:

Δu=∇2u=∇⋅∇u=0=1r2∂∂r(r2∂u∂r)+1r2sin⁡θ∂∂θ(sin⁡θ∂u∂θ)+1r2sin2⁡θ∂2u∂2φ

对其分离变量 u(r,θ,φ)=R(r)Y(θ,φ), 容易得到

1Rddr(r2dRdr)=−1Y(1sin⁡θ∂∂θ(sin⁡θ∂Y∂θ)+1sin2⁡θ∂2Y∂ϕ2)=λ

2 径向衰减

常常选取 λ=l(l+1), 于是

ddr(r2dRdr)−l(l+1)R=0

这是欧拉方程, 根据换元 r=es 得到它的通解是

R(r)=Arl+Br−(l+1),l≥0

通常情况下, 场不会传播到无穷远处, 也就是在解 R(r) 的定义域包含 r=∞ 时因为 rl→∞ 导致 R(r)→∞ 而不满足物理意义, 所以通常取 A=0, 这将只剩下 r−l−1 的项, 表明场随着 r 增加而逐渐衰减. 同理, 在解的定义域包含 r=0 时, 通常取 B=0.

3 球面谐波

3.1 球谐函数的引出

由上述已知, 拉普拉斯方程的球面偏导部分是

1sin⁡θ∂∂θ(sin⁡θ∂Y∂θ)+1sin2⁡θ∂2Y∂ϕ2+l(l+1)Y=0

进一步分离变量 Y(θ,φ)=Θ(θ)Φ(φ), 有

(1){d2Φd2φ+m2Φ=0Φ(φ)=Φ(φ+2π)(2)sin⁡θddθ(sin⁡θdΘdθ)+[(l+1)lsin2⁡θ−m2]Θ=0

方程(1)没什么好说的, 解是 Φm(φ)=eimφ. 对方程(2)做换元 θ=arccos⁡x, 得到

(3)(1−x2)d2Θdx2−2xdΘdx+[l(l+1)−m21−x2]Θ=0

这被称为连带勒让德方程. 特别的, 当 m=0 时, 上式退化为勒让德方程:

(4)Θ→Pl:(1−x2)d2Pldx2−2xdPldx+l(l+1)Pl=0

该方程的解就是勒让德函数 Pl(x).

对连带勒让德方程(3)做换元 y(x)=Θ(x)(1−x2)−m2, 可得到

Θ=(1−x2)m2yΘ′=(1−x2)m2y′−m(1−x2)m2−1xyΘ″=(1−x2)m2y′′−2m(1−x2)m2−1xy′−m(1−x2)m2−1y+m(m−2)(1−x2)m2−2x2y(5)⇒(1−x2)y′′−2(m+1)xy′+[l(l+1)−m(m+1)]y=0

将该方程与勒让德方程(4)做比较, 先对勒让德方程求导 m 次, 记 Pl(m) 表示对 Pl 求 m 次导, 利用莱布尼兹公式 (uv)(m)=∑k=0m(mk)u(m−k)v(k), 可以得到

0=[(1−x2)Pl(2)−2xPl(1)+l(l+1)Pl](m)=0=[(1−x2)Pl(m+2)−2mxPl(m+1)−m(m−1)Pl(m)]−2[xPl(m+1)+mPl(m)]+l(l+1)Pl(m)=(1−x2)Pl(m+2)−2(m+1)xPl(m+1)+[l(l+1)−m(m+1)]Pl(m)

这恰好和式(5)的形式一致, 也就是说方程(5)的解是

ylm=dmdxmPl(x),0≤m≤l

那么连带勒让德方程(3)的解记为连带勒让德函数 Plm(x), 回忆前面做过换元 y(x)=Θ(x)(1−x2)−m2, 那么

Plm(x)=(1−x2)m2dmdxmPl(x),0≤m≤l

所以, 最后的球谐函数就是

Ylm(θ,φ)=Plm(cos⁡θ)Φm(φ)

如果要归一化 Y, 也就是令 ∫02πdφ∫0πYlm(θ,φ)sin⁡θdθ=1, 则可以得到归一化的球谐函数:

Ylm(θ,φ)=12π2l+12(l−m)!(l+m)! Plm(cos⁡θ)eimφ

同时也能够给出 −l≤m<0 的项满足,

Pl−m(x)=(l−m)!(l+m)!Plm(x)Yl−m(θ,φ)=Ylm∗(θ,φ)

于是 ∇2u=0 的通解表示为

(6)u(r,θ,φ)=R(r)Ylm(θ,φ)=∑l=0∞∑m=−ll(Arl+Br−(l+1))12π2l+12(l−m)!(l+m)! Plm(cos⁡θ)eimφ

另外, 这里列出连带勒让德函数的前几项:

P00(x)=1,P10(x)=x,P11(x)=1−x2,P20(x)=12(3x2−1),P21(x)=3x1−x2,P22(x)=3(1−x2),P30(x)=12(5x3−3x),P31(x)=32(5x2−1)1−x2,P32(x)=15x(1−x2),P33(x)=15(1−x2)3/2.

在量子力学中,经常在连带勒让德函数前面加上系数 (−1)m 变成:

Plm(x)=(−1)m(1−x2)m2dmdxmPl(x),0≤m≤l

你会发现这个新的系数实际上不影响连带函数的正交完备性,而完全只是量子力学的一个相位规定。

3.2 Y1m 和偶极场

在静磁学中, 如果场点不存在电流(∇×B→=0), 那么磁感应强度 B→ 可以表示成磁标势 Φm 的梯度 B→=−∇Φm, 并且由于 ∇⋅B→=0, 所以磁标势自动满足拉普拉斯方程 ∇2Φm=0. 在其通解形式中, 前几项球谐函数有显著的物理意义.

上述的归一化球谐函数中, Y10(θ)=34πcos⁡θ∝cos⁡θ=zr, R1(r)=Ar+Br2=Br2 (取 A=0). 于是 Φm 的一个解可以写作

Φpole(r,θ,φ)=Cr2cos⁡θ,C=const

求梯度得到

B→=−(∂rΦpole, 1r∂θΦpole, 0)=C⋅(2cos⁡θr3, sin⁡θr3, 0)

回忆磁偶极子在远处的磁场表达式:

B→=μ04π3(n→⋅m→)n→−m→r3=μ04π3cos⁡θe→r−e→zr3=μ04π3cos⁡θe→r−(cos⁡θe→r−sin⁡θe→θ)r3=μ04π2cos⁡θe→r+sin⁡θe→θr3

发现两式形式一致, 只要取 C=μ04π, 就可以用球谐函数 Y10 描述偶极磁场.

另外, 由于 Y1±1 是虚数, 所以为保证物理意义, 选取它们的实数线性组合 Y1x=Y11+Y1−12 和 Y1y=Y11−Y1−1i2, 显然 Y1x∝sin⁡θcos⁡φ=xr, Y1y∝sin⁡θsin⁡φ=yr 为实数. 对于三维旋转, 考察球谐函数的变换, 比如在直角坐标系中, 绕 y 轴旋转的变换矩阵是

Ry(β)=(cos⁡β0sin⁡β010−sin⁡β0cos⁡β)

变换导致 (x,y,z)→(x′,y′,z′),(r,θ,φ)→(r′,θ′,φ′), 即:

θ′=arccos⁡(z′r)=arccos⁡(−xsin⁡β+zcos⁡βr)=arccos⁡(−sin⁡θcos⁡φsin⁡β+cos⁡θcos⁡β)φ′=arctan⁡y′x′=arctan⁡yxcos⁡β+zsin⁡β=arctan⁡sin⁡θsin⁡φsin⁡θcos⁡φcos⁡β+cos⁡θsin⁡β

如果取 β=π2, 就有

θ′=arccos⁡(−sin⁡θcos⁡φ)φ′=arctan⁡(sin⁡θsin⁡φcos⁡θ)

也就是说, 在绕 y 轴旋转 90∘ 后, Y1x 变成了描述朝向 z 轴的偶极子 Y10 :

Y1x(θ′,φ′)∝sin⁡θ′cos⁡φ′=sin⁡(arccos⁡(−sin⁡θcos⁡φ)) ⋅ cos⁡(arctan⁡(sin⁡θsin⁡φcos⁡θ))=1−(sin⁡θcos⁡φ)2  ⋅11+(sin⁡θsin⁡φcos⁡θ)2=1−(sin⁡θcos⁡φ)2  ⋅cos⁡θcos2⁡θ+(sin⁡θsin⁡φ)2=cos⁡θ∝Y10(θ)

因此, Y1x 实际上描述的是沿 x 轴方向的磁偶极子, 类似的知道 Y1y 描述沿 y 轴方向的磁偶极子. 在量子力学里面, 上面的旋转操作可以用 Winger−D 矩阵作用于球谐函数得到, 这个矩阵的元素很丑陋, 可以在 MMA 中直接调用 Winger−D 函数快速计算.

3.3 Y2m 和四极场

同理, 对于归一化球谐函数 Y20(θ)=516π(3cos2⁡θ−1), 磁标势 Φquad(r,θ)=R⋅g20(Rr)3516π(3cos2⁡θ−1), g20 为常数, 计算得到

Br=−∂Φ∂r=3g20(R4r4)516π(3cos2⁡θ−1)Bθ=−1r∂Φ∂θ=3g20(R4r4)516πsin⁡2θBφ=0

关于球谐系数 glm 的具体拟合值, 可在第 14 版国际参考地磁场(IGRF)的官方文件 中查阅. 对于 m≠0 的球谐函数, 同样选取他们的实数线性组合, 可以得到他们表示的四极场只是 Y20(θ) 的旋转, 但是形式更加复杂.

4 球谐系数

4.1 函数 1|r→−a→| 的勒让德系数

处理电磁场的问题中常常遇到函数 f(r)=1|r−a|,出于简洁,首先假设 r=(r,θ,φ),a=(a,0,0),那么可见函数 f(r) 是轴对称的(即与 φ 无关),可展开为:

1|r−a|=∑l=0∞Al(r,a)Pl(cos⁡θ)

当 θ=0:

1|r−a|=1|r−a|={1r11−(a/r)=1r∑l=0∞(ar)l=∑l=0∞alrl+1,r>a1a11−(r/a)=1a∑l=0∞(ra)l=∑l=0∞rlal+1,r<a

简记为

1|r−a|=∑l=0∞r<lr>l+1,r<=min(r,a), r>=max(r,a)

由 Pl(cos⁡0)=Pl(1)=1 和对比系数:

1|r−a|=∑l=0∞Al(r,a)Pl(cos⁡0)=∑l=0∞r<lr>l+1

可知 Al(r,a)Pl(cos⁡0)=r<lr>l+1:

1|r−a|=∑l=0∞r<lr>l+1Pl(cos⁡θ)

现在旋转得到 a=(a,θ′,φ′),由于 |r−a| 只取决于两个矢量之间的相对夹角 γ=arccos⁡r⋅a|r⋅a|,所以对于一般情形,有

1|r−a|=∑l=0∞r<lr>l+1Pl(cos⁡γ),γ=arccos⁡r⋅a|r⋅a|, r<=min(r,a), r>=max(r,a)

或者,也可以从勒让德多项式的母函数表示来证明上述结果。

4.2 已知球对称边界条件

现在需要确定系数 Blm, 假设 A=0, 通解写成:

u(r,θ,φ)=∑l=0∞∑m=−llBlmrl+1Ylm(θ,φ)

其中 Ylm(θ,φ)=(2l+1)4π(l−m)!(l+m)!Plm(cos⁡θ)eimφ, 它是正交归一化的:

∫02π∫0πYlm(θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ=δll′δmm′

将通解两边乘以 Yl′m′∗(θ,φ), 然后积分:

∫02π∫0πu(a,θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ=∑l=0∞∑m=−llBlmal+1∫02π∫0πYlm(θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ

得到上式需要假定边界条件为 r=a 时 u(r,θ,φ)=u(a,θ,φ) 已知。利用正交归一关系, 右边只剩下 l=l′ 且 m=m′ 时的项:

∫02π∫0πu(a,θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ=Bl′m′al′+1

从而可解出:

Bl′m′=al′+1∫02π∫0πu(a,θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ

同理,假设 B=0,可以得到:

Al′m′=a−l′∫02π∫0πu(a,θ,φ)Yl′m′∗(θ,φ)sin⁡θdθdφ

4.2.1 反解偶极场

举个例子,假设已知边界电势 φ(a,θ,φ)=V0cos⁡θ,求解这是一个什么场。 由于分布与 φ 无关,只有 m=0 的项贡献。根据球谐函数的定义:

Y00=14π,Y10=34πcos⁡θ

由此可见,我们的边界分布可以写为:

φ(a,θ)=V04π3Y10(θ,φ)

所以球谐系数根据前面提到的公式计算:

Blm=al+1⋅V04π3∫ΩY10Ylm∗dΩ

得到 B10=V0a24π3,即 φ(r,θ)=V0a2r2cos⁡θ  (r≥a),这恰好是电偶极矩为 p=4πϵ0V0a2 的点电偶极子的电势(或者前面提过的磁偶极子)。